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ABSTRACT 

We use a particle tagging technique to dynamically populate the ./V-body Via Lactea II high- 
resolution simulation with stars. The method is calibrated using the observed luminosity function 
of Milky Way satellites and the concentration of their stellar populations, and self-consistently follows 
the accretion and disruption of progenitor dwarfs and the build-up of the stellar halo in a cosmo- 
logical "live host" . Simple prescriptions for assigning stellar populations to collisionless particles are 
able to reproduce many properties of the observed Milky Way halo and its surviving dwarf satellites, 
like velocity dispersions, sizes, brightness profiles, metallicities, and spatial distribution. Our model 
predicts the existence of approximately 1,850 subhalos harboring "extremely faint" satellites (with 
mass-to-light ratios > 5 x 10 3 ) lying beyond the Sloan Digital Sky Survey detection threshold. Of 
these, about 20 are "first galaxies", i.e. satellites that formed a stellar mass above 10 M Q before red- 
shift 9. The ten most luminous satellites (L > 10 6 Lq) in the simulation are hosted by subhalos with 
peak circular velocities today in the range V max = 10 — 40 kms" 1 that have shed between 80% and 
99% of their dark mass after being accreted at redshifts 1.7 < z < 4.6. The satellite maximum circular 
velocity Vmax and stellar line-of-sight velocity dispersion <7i os today follow the relation F max = 2.2<ti os . 
We apply a standard mass estimation algorithm based on Jeans modelling of the line-of-sight velocity 
dispersion profiles to the simulated dwarf spheroidals, and test the accuracy of this technique. The 
inner (within 300 pc) mass-luminosity relation for currently detectable satellites is nearly flat in our 
model, in qualitative agreement with the "common mass scale" found in Milky Way dwarfs. We do, 
however, predict a weak, but significant positive correlation for these objects: M300 oc £°- 088±0 - 024 . 
Subject headings: Galaxy: halo - kinematics and dynamics - galaxies: dwarf - formation - structure 
- methods: numerical 



^1 _ 1. INTRODUCTION 

£> ' Galaxy halos are one of the crucial testing grounds 
C*~) \ for ACDM structure formation scenarios, as they con- 
OO . tain the imprints of their assembly via the hierarchi- 
' cal merging and accretion of many smaller p rogeni- 
*0 I tors (jWhite fc ReesHl978t iBhimenthal et al.lll984l) . Sub- 
■ units are predicted to collapse at high redshift and have 
C^) ' high density concentrations, and when they merge into 
more massive hosts they are able to resist tidal disrup- 
tions and survive in large numb ers as bound substruc- 
tures (e.g. iDiemand et al.ll2008 l). The dwarf spheroidal 
(dSph) satellites present today in the halos of massive 
, galaxies like our own Milky Way (MW) are thought 
5h ' to represent a fossil record of such early merging ac- 
tivity. The old tally of MW "classical" dwarfs with 
Mv ^ ~ 10 (e.g. lMatedll998l ) has been recently revised 
by the discovery of a new, large population of ultra- 
faint dwarfs in the Sloan Digital Sky Surv ey (SDSS) 
(jWillman et al.l 12001 iBelokurov etail [2007h . Follow- 
up spectroscopic observations have shown that these 
extreme objects are the most dark matter-dominated 
galaxies as well as some of the most metal-poor stel- 
lar systems known in the universe jSimon fc Oehall2007t 
iKirbv et alj|2008t iGeha et al.l 120091 ) . While most of the 
observed dSphs seem to be in pressure-supported equilib- 
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rium , others (like the Sagitt arius dwarf, see llbata et al.1 
120011; iMaiewski et al.l 120031) are being ripped by tidal 
forces into long streams of stars that will eventually 
dissolve into the Galactic stellar halo. Large coher- 
ent st ellar streams have been ob served in Andromeda's 
halo (McCon nachie et al.1 l2009[) and in the outskirts 



of other nearby galax ies (jMartmez-Delgado et al.l [20101 : 
iMouhcine eTalBoToT) . 

The predicted counts of bound substructures vastly 
exceed the number of known satellites of the MW, 
a "missing satellite problem" that has been the sub- 
ject of many numerical studies over the last decade 
Ce.g.lKlvpin et al.lll999t iMoore et aLlll999t IMadau et al.1 
2008a|P While this discrepancy is partially alleviated 
by correcting the data to account for the incomplete 
sky coverage and distance-magnit ude or surface bright- 
ness limits of curr ent surveys (IKoposov et al.1 120081 : 
iTollerud et all [20081 : iBullock et al.ll2010D . it is~generallv 
agreed that a complete solution to the apparent con- 
flict between the Galaxy's relatively few observed stel- 
lar clumps and the abundance of halo substructure in 
cold dark matter simulations requires either a modifi- 
cation of the nature of the dark matter p article that 
would act to suppress small-scale power (e.g. lColm et al.1 
120001) or an extremely low star formation efficiency in 



small dark matter subunits fe.g. iKauffmann et al.[[l993 
Benson et alj|2002t iBTuTock et al.ll2000t: lSomervillell2002 



Kravtsov et al.l l2004t IMadau et al.ll2008bt iKuhlen et al.1 
20111) . 



Measurements of the line-of-sight velocity dispersion 
as a function of radius in dwarf satellites, combined 
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with a spherical Jeans analysis, may have provided an- 
other "small-scale puzzle". Despite spanning about 5 
decades in luminosity, MW dSphs appear to contain a 
similar amount of dark mass, M 300 ~ 10 7 Mq, within 
a fixed aperture radius of 300 pc from their centers 
(jStrigari et all 120081 : iWalker et aLl I2009Q . It is unclear 
at this stage whether the flat inner mass-luminosity rela- 
tion inferred, M300 oc 2 y 03±0 03 j is an accident of small 
number statistics, an artifact of the spherical Jeans anal- 
ysis technique, arises from an observational selection bias 
([Bullock et al.l[20To| ). is ev idence for a new fu ndamental 
scale in galaxy formation (jMunoz et al.l[2009[). or arises 
naturally in the ACDM paradigm (e.g. iMaccio et aTl 
[2001 IStringer et alJl20Toh . 

Numerical simulations of the assembly of the halo of 
Milky Way-sized galaxies have the potential to shed light 
on the nature of dSph galaxies and test the ACDM hi- 
erarchical clustering paradigm. Owing to their compu- 
tational expense, state-of-the-art cosmological hydrody- 
namic simulations that include gas cooling, star forma- 
tion, and supernova feedback, and can self-consistently 
follow the formation and evolution of MW substructure 
down to the present e poch, are limited to dark particle 
masses m p > 10 5 M Q (jGuedes et al.ll20ilT ). and can well 
resolve only the most massive subhalos. To be able to 
reach the necessary resolution, we use in this paper a 
particle tagging technique to dynamically populate the 
A^-body Via Lactea II (hereafter VLII) high-resolution 
simulation with stars, and self-consistently follow the 
accretion and disruption of progenitor dwarfs and the 
build-up of the stellar halo in a cosmological "live host" . 
We show how simple prescriptions for assigning stellar 
populations to collisionless particles are able to repro- 
duce many properties of the observed MW halo and its 
surviving dwarf satellites, like velocity dispersions, sizes, 
brightness profiles, metallicities, and spatial distribution. 
More importantly, with over 1.6 million stellar particles 
tagged in post-processing, we can then mimic the ob- 
servations of the line-of-sight velocity dispersion profiles 
of the simulated dSphs, apply a standard mass estima- 
tion algorithm based on spherical Jeans modelling, test 
the accuracy of this technique, and confirm or dispel the 
existence of a common mass scale for satellite galaxies. 

A number of previous implementations of the 
particle-tagg ing approach e x ists i n the litera- 
turc (e.g. IWhite fc Springell WM. iDiemand et all 



TABLE 1 
Build-up of the VLII stellar halo 



20051: iBullock fc Johnston! 120051:" 



Cooper et al.l 120101 : lTumlinsonll201oi ). The approxima 



Moore et all 120061 : 



tion breaks down when baryons become dynamically 
important, i.e. in the inner regions of the MW and in 
the galactic disk. The orbits of some subhalos will pass 
through these regions, and one has to keep in mind and 
eventually correct for the fact that some substructure 
will actually feel the potential of baryonic configurations. 
Dwarf spheroidal satellites on the other hand seem to be 
dark matter-dominated on all scales, and our simulated 
subhalo potentials can be directly compared to the 
observed stellar kinematics. VLIFs high mass resolution 
allows the survival of small dark matter subhalos long 
after they are accreted by the main host. Resolved at 
the present epoch with anywhere between 1,000 and 
100,000 particles, these systems, in combination with 
an appropriate technique to designate stellar particles, 
serve as a testbed for the techniques commonly used to 
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Note. — Columns 1 and 2 give the VLII snapshot number and 
its corresponding redshift. Columns 3 and 4 list the number of 
subhalos that reached their maximum dark matter mass in 
that snapshot, together with the total stellar mass they contribute 
to the VLII main host at z = 0. 

estimate the masses of dSph galaxies from real data. 

The paper is organized as follows. In § 2 we describe 
the VLII simulation and particle tagging technique, to- 
gether with the resulting stellar halo and the population 
of dwarf satellites surviving today. In § 3 we describe 
the Jeans analysis used to estimate masses in pressure- 
supported collisionless systems, and assess its perfor- 
mance in correctly estimating the underlying mass distri- 
butions of the simulated dSph galaxies. We summarize 
our findings in § 4. 

2. SIMULATION AND PARTICLE TAGGING TECHNIQUE 

The VLII simulation, one of the highest-precision cal- 
culations of the assembly of a Galact ic halo to date 
(|Diemand et al.ll2008l iZemp et al.|l2009fl was p erformed 
with the PKDGRAV tree-code istadej |200l| . Initial 
conditions were based on the WM AP 3-year data re- 
lease cosmology (jSpergel et al.ll200"7h and were generated 
with a modified, paralle l version of the GRAFIC2 pack- 
age (jBertschingerl l200ll ). VLII employs just over one 
billion 4, 100 M Q particles to model the formation of a 
M 2 oo = 1-9 x 10 12 M Q Milky- Way size halo and its sub- 
structure. About 40,000 subhalos of masses as low as 
10 6 Mq are resolved today within the host's r2oo = 402 
kpc (the radius enclosing an average density 200 times 
the mean matter value). 

The particle tagging technique makes full use of 27 
snapshots (out of 400 total) of the simulation, spanning 
the assembly history of the host between z = 27.54 and 
the present. Numbered from to 26, the snapshots are 
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Fig. 1. — Distributions of the peak mass at infall for 3,204 VLII subhalos (in blue) above = 10 7 Mq (left panel) and of their infall 
redshift (right panel). Tagged subsets of dark matter particles in these systems are followed in post-processing from infall to the present. 
Subhalos that at z = host a SDSS-detectable dSph populate the higher-mass and lower-redshift ends of these distrubutions (in red). The 
gaps in the redshift distribution correspond to Az = 1 intervals where VLII snapshots were not analyzed in detail (see Table [TJ. 



listed in Table[JJtogether with the corresponding redshift. 
In each snapshot, all ( sub)halos were identi fied with the 
6DFOF group finder (jDiemand et all 120061) . and linked 
from snapshot to snapshot to their most massive progen- 
itor: the halo 'tracks' built in this way contain all the 
time-dependent structural information necessary for our 
study. We use them to identify the simulation snapshot 
in which each subhalo reaches its maximum mass Mh be- 
fore being accreted by the main host and stripped. By 
neglecting all subhalos with Mh < 10 7 M Q , we restrict 
our analysis to 3,204 such tracks that contribute mass to 
the main host halo at z = 0. The mass distribution of the 
selected subhalos at infall and the distribution of infall 
redshifts are shown in Figure [TJ The gaps in the red- 
shift histogram correspond to epochs in which no VLII 
snapshots were studied in detail (see Tabic [JJ . Due to 
the great computational effort necessary in the analysis 
of a single VLII snapshot, 6DFOF was applied only to 
the 27 ones listed here, and expanding this analysis to 
more snapshots is beyond the scope of our work. Ha- 
los that are accreted by the main host in these intervals 
are tagged at the next available snapshot. As the halos 
we consider are composed of a large number of particles, 
they are well resolved even ~0.5-0.8 Gyr (typical snap- 
shot separation) after accretion, making it unlikely that 
this time discretization affects significantly our results. 

A subset of dark matter particles in each subhalo is 
then "tagged" as stars at infall, and no subsequent star 
formation is assumed to take place - hydrodynamical 
simulations suggest that ram pressure stripping would 
deplete the subhalo gas reserv oir and shut star forma- 
tion off shortly thereafter (e.g. iMaver et al.l[20 06). It is 
therefore at infall that tagged particles acquire a stellar 
mass m sp , an age i sp , and a mctallicity [Fc/H] sp : any 
subsequent evolution of the tagged stellar populations 
is purely photometric and kinematical in character, as 
systems age and are accreted and disrupted in a cosmo- 
logical "live host" . 

The following three prescriptions determine m sp , t sp , 
and [Fe/H] S p, and were tuned to match certain structural 
and statistical properties of the surviving dSph popula- 



tion at z = 0: the star formation efficiency M^/Mh is 
fit to reproduce the satellite luminosity function, the se- 
lection of particles painted in halos is made to fit the 
current extent of their luminous component, and the as- 
signed metallicities are designed to place satellites on the 
observed metallicity-luminosity relation today. The first 
two prescriptions thus determine the value m sp for each 
tagged particle, while the third determines its [Fe/H] sp . 
The age t sp of the stellar population consistently follows 
from the star formation efficiency prescription and the 
history of dark matter build-up of the subhalo it was 
tagged in. 

In t he spirit of IKopqsov et al.l (|2009f ) and IKravtsovj 
(|2010[) . the total stellar mass tagged at infall is assumed 
to scale as a power-law of the host subhalo mass, 

M h \10 9 M Q J 1 ' 

IKravtsovj (|2010D used a slope of 1.5 but argued that 
the results would not change drastically if a somewhat 
steepe r (2.0-2.5) slope was assumed. iKoposov et al.l 
(|2009t) similarly adopted a slope of 2.0. Our relation has 
a somewhat higher normalization than both. The differ- 
ence likely arises from the conversi on of the assign e d stel - 
lar mass to a V-band luminosity. IKoposov et al.l ()2009l ) 
assume a mass-to-light ratio M/L = 1 to be typical for 
old stellar populations, while wc find that M/L ~ 3 — 4 
better fits the range of ages and metallicities we assign to 
the infalling s atellites (see bel ow). A similar argument 
holds for the iKravtsovl (|2010l ) relation, which directly 
assigns a V-band luminosity to a halo of a given mass, 
instead of a stellar mass. Additionally, some of the halos 
in VLII lose a fraction of their stellar mass after accre- 
tion from tidal stripping, an effect that may not have 
been treated in the other two studies. 

Our stellar mass-subhalo mass relation at infall has a 
scaling similar to the low-mass end (Mh ~ 10 11 ' 5 M Q ) 
galaxy stellar mass-halo virial mass relation, M^/Mh oc 
M^ 9 , derived at z = 1 by iBehroozi et all (|2010l ) using 
an abundance matching technique. It implies a star for- 
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mation efficiency (the fraction of baryons turned into 
stars), /* = (nM/Qb)(M*/Mh), that is close to 1% in 
Mh = 10 10 M0 subhalos, and becomes negligibly small 
in systems below Mh = 10 s M . The detailed physical 
explanation of this behavior is still an area of active re- 
search. Here we stress that it is the simple prescription 
in equation (JT|) that largely determines the luminosity 
function of surviving dwarf satellites today and empiri- 
cally bypasses the missing satellite problem by strongly 
suppressing star formation in low-mass systems. 

The total stellar mass within an infalling subhalo is 
equally distributed to its 1% most tightly bound dark 
matter particles. This number determines the concen- 
tration of the stellar system at the time of infall, and 
therefore governs the amount of stellar material stripped 
at later times as well as the present-day structural prop- 
erties of surviving satellites, li ke their half-li g ht rad ii and 
central surface brightnesses. iCooper et al.l (|2010f) have 
shown that tagging the 1% highest total binding energy 
particles with stellar populations provides a good fit to 
the distribution of half-light radii in MW dSphs. Our 
analysis below confirms this result. 

The metallicity assigned to each particle at infall is also 
a function of Mh, 

[Fe/H] sp = -7.87 + 0.9 x log (^^) , (2) 

i.e. more massive subhalos are assumed to retain more 
of their enriched material and have a higher [Fe/H] than 
smaller systems. All particles initially tagged in the same 
subhalo are assigned the same metallicity. This prescrip- 
tion does not account for the timescale of enrichment 
of the stellar systems with iron or alpha-elements. In 
the present work metallicity only enters in the calcula- 
tion of the luminosity of a given stellar system, and docs 
not affect our kinematic studies. Two processes turn the 
above scaling into th e observed present -day luminosity- 
metallicity relation of lKirbv et al.l (|2008l ): (i) stellar mass 
loss from tidal stripping and (ii) dimming of the stel- 
lar populations with age (see the right panel of Fig. |3]). 
When comparing the above relation to the one derived for 
SDSS gala xies, we note that it li es below the extrapola- 
tion of the lTremonti et a l. (2004) relation. This suggests 
that processes governing metal formation and retention 
in dSphs may be significantly different from those in more 
massive galaxies. 

Finally, while all particles in a given subhalo are tagged 
at a single epoch (infall) , they have a distribution of stel- 
lar ages resulting from applying equation ([T]) to the his- 
tory of dark matter build-up of each halo prior to its infall 
and assigning the respective age t sp to the corresponding 
fraction of tagged stellar mass. This process mimics the 
extended star for mation histories of dwarf satellites in 
the MW halo fsee lOrban etTaTl [20081) . 

After all stellar particles are tagged with the above 
properties, they are traced down to the z = snapshot 
in the simulation. Some arc still part of surviving dSph 
satellites; many form coherent stellar streams and con- 
tribute to the smooth stellar halo component. Their fi- 
nal V-band luminosities arc sclf-consistcntly determined 
for their age, stellar mass, and metal licity from the stel- 
lar po pulation synthesis models of iBruzual &: Chariot! 
(|2003l) (assuming a Chabrier initial mass function, 



IChabriedl2003l) . 

2.1. Assembly history 

A total of 1.6 million tagged particles end up in the 
main host at redshift z = (0.35% of all particles within 
r 2 oo): most are tidally stripped from their subhalos af- 
ter infall, and only 52% of them remain bound to the 
1,925 self-gravitating clumps that partially survive the 
accretion event. The total mass of the stellar halo to- 
day is M* = 1.6 x 10 9 Mq within 402 kpc (of which 
2.3 x 10 s M Q is in bound substructures), and Table Q] 
lists the contribution to from each snapshot: the 
peak contribution comes from redshift 2. The mass in 
the stripped component between 1 kpc and 40 kpc is 
5 x 10 s M Q , comparable to recent estimates for th e MW 
halo mass within the same range dBell et al.ll2008l ). 

In terms of its mass assembly history, the stellar com- 
ponent of the VLII halo is already in place by redshift z = 
1.5, while the dark component continues to grow from 
mergers with smaller subunits down to the present epoch. 
The faster assembly of the stellar component follows from 
the very steep mass-star formation efficiency relation 
adopted in our model. A similar behavior was seen in 
the assembly history of the six Aquarius stellar halos by 
ICooper et all (|2010t i. About 50% of today's stellar halo 
can be traced back to a subhalo of M* = 8.8 x 10 8 M© 
that fell in at z = 2. The three next most massive con- 
tributors, with M*/ M Q = (4.2 x 10 8 , 1.2 x 10 8 , 4.5 x 10 7 ), 
were accreted at redshifts z =7.8, 3.2, and 1.6, respec- 
tively. Combined, these four progenitor systems account 
for about 88% of the total stellar halo. 

Figure [2] shows an image of the projected VLII dark 
matter and stellar halo mass density at 3 different red- 
shifts. The disruption of the stellar component of an 
individual satellite is depicted in the right panels, where 
particles initially tagged in one subhalo at infall (the sys- 
tem reaches its Mh at z = 2) are violently stripped by 
tides and form shell- and stream-like structures. The 
bottom left and middle panels make it clear how the 
adopted steep stellar mass-halo mass relation offers an 
empirical solution to the "missing satellite problem" . 
While the dark halo exhibits a vast amount of surviving 
self-bound substructure, only the most massive subhalos 
form enough stellar mass at infall to still be visible in the 
stellar halo today. 

2.2. Surviving satellites 

Surviving satellites in the z = VLII snapshot were 
identified by considering all particles that fall within the 
tidal radius of a bound subclump (identified in the dark 
matter density field by the 6DFOF group finder, see 
IDiemand et al.ll2006D . We then performed a peculiar ve- 
locity cut of V p < 2V max in order to remove particles 
belonging to the main halo that just happen to be flying 
through the given satellite at that time. The total lu- 
minosity of each satellite was computed by summing up 
the contributions of all stellar particles linked to it, and 
its metallicity determined from their mass- weighted aver- 
age. Satellites of absolute magnitude M v located closer 
than the SDSS completeness limit r max (jTollerud et al.1 
2008) fr om the main halo center, 

W = (—1—) V3 X W (-™^- 5 .23)/3 MpCj (3) 

\47r/ DR5 / 
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Fig. 2. — Projected mass surface density of the VLII dark matter halo (left panels) and its stellar halo (central panels), at redshift 2 (top), 
1 (middle) and today (bottom). The image size is 400 physical kpc in all panels. For comparison, the T20Q radius of the VLII dark matter 
halo is 94 kpc, 180 kpc, and 400 kpc at these 3 epochs, respectively. The color scale (from dark to bright) ranges from 2 X 10 s Mq kpc -2 
to 10 9 Mq kpc -2 for the dark matter, and from 10 — 2 Mq kpc -2 to 2 X 10 8 Mq kpc~ 2 for the stars. Right panels: the disruption of the 
stellar component of an individual satellite tagged at infall (z = 2). 
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were flagged as observable in current surveys (/dr.5 = 
0.194 here is the sky covering fraction of the SDSS DR5). 
Only 65 of the 1,925 surviving satellites present today in 
the VLII stellar halo are bright enough to be detected by 
an all-sky SDSS-like survey, out to a galactocentric dis- 
tance o f 280 kpc, corres ponding to the Milky Way virial 
radius (|Xue et all 120081) . Note that, as the SDSS detec- 
tion efficiency is nearly independent of central surface 
brigh tness for fj, y 30 mag arcsec -2 (|Koposov et al.1 
120081 ). no surface brightness threshold has been applied 
to our sample. We have checked that all but 13 of the 
observable VLII satellites have a central surface bright- 
ness of at least 30 mag arcsec -2 . The 13 low surface 
brightness systems are all faint (My > —5) and have 
half-light radii r c g close to the gravitational softening 
scale of the simulation, i.e. they have been artificially 
puffed up. We have estimated that, without such a nu- 
merical effect, these subhalos would actually lie above 
the SDSS surface brightness limit. Our prescriptions do 
not therefore produce a large population of low surface 
brightness "stealth" satell ites awaiting d i scove ry in the 
MW halo, as sugg ested bv lBullock etall pOlOh . 

The stellar masses today of the observable VLII satel- 
lites range from 2,000 M Q to 10 8 M , and their dark 
masses from 4 x 10 6 M Q to 3 x 10 9 M . They have mass- 
to- light ratios that range from about 20 to 10 5 , and infall 
redshifts mostly between 1.5 and 8 (with a few falling in 
later, between redshifts 0.6 and 1.1). Their mass-to-light 
ratios have decreased by a factor ranging between a few 
and 20, as their extended dark matter halos have been 
tidally stripped more than their compact luminous com- 
ponents. The ten most luminous satellites (L > 10 6 Lq) 
in the simulation are hosted by subhalos with peak cir- 
cular velocities today in the range V max = 10 — 40 kms -1 
that have shed between 80% and 99% of their dark mass 
after being accreted at redshifts 1.7 < z < 4.6. Their 
Knax values have decreased by a factor of 2.5 on the av- 
erage (and up to 6 in individual cases) after accretion, 
so their present-day Kn ax is not a good indicator of the 
subhalo maximum mass. 

The differential satellite luminosity function is shown 
in the left panel of Figure [3l plotted against the lu- 
minosity function of observed MW dwarfs (see also 
w orks bvlKoposovet al l 120091: iMa ccio et al.l l2010f). As 
in lKoposov et al.l ~ 2009(1 . we use different scales for the 
number of satellites brighter and fainter than My = — 11, 
introduced to produce a continuous luminosity function 
despite the limited SDSS DR5 sky coverage of /dr,5 = 
0.194. The shaded region to the left of My = — 11 spans 
the la extent in the number of faint, but detectable satel- 
lites that fall within the SDSS DR5 footprint in each of 
2,000 randomly oriented sky projections of the VLII stel- 
lar halo. In the case that the satellites are isotropically 
distributed, our exercise should produce a result equiva- 
lent to a scaling of the number counts of observed faint 
satellites by the sky coverage fraction of the SDSS survey. 
Satellites to the right of the dashed line, corresponding 
to My = — 11 are always detectable. The shaded region 
to the right of My = — 11 covers the ltr Poisson errors in 
each bin there. The lack of Magellanic Cloud equivalents 
in VLII is responsible for the discrepant bright end of 
the luminosity function, as the simulation does not con- 
tain satellites more luminous than My < —14. Studies 
suggest that satellites this bright are rare around Milky 



Way- like galaxies (|Bovlan-Kolchin et al.|[20lol : iLiu et ail 
120111) . In the right panel of the same figure, we plot 
the luminosity-metallicity relatio n of surviving sate llites 
compared to data for MW dSphs (jKirbv et al.ll2008l and 
references therein) and their best-fit model to the data. 
The relation assumed at infall (eq. [2]) has been modified 
by tidal mass losses and by the aging of stellar popula- 
tions. 

In FigurelH the cumulative galactocentric distance dis- 
tribution of VLII luminous subhalos is compared to that 
of MW satellites fainter than My > — 12. Even though 
the final spatial distribution of dwarfs depends on their 
accretion and orbital history, the two distributions track 
each other reasonably well, though VLII shows an excess 
of close-in satellites compared to the observations. The 
figure also shows (right top panel) the present-day stellar 
mass-subhalo mass relation. Much like the mctallicity- 
luminosity relation in Figurc[3l it has been modified from 
its infall value by tidal mass losses. As accreted sys- 
tems lose a disproportionate amount of their dark mass, 
the apparent (inferred today) star formation efficiency of 
dSphs may exceed the true one by as much as a factor 
of 100. We also explore (bottom left panel of Fig. |4j 
how the average 1-D line-of-sight velocity dispersion a\ os 
of dSphs (an observable quantity) compares to the peak 
(3-D) circular velocity F max of the system's host halo (a 
quantity of theoretical interest). Early attempts at this 
conversion assumed a proportionality facto r of y/3 mo- 
tivate d by the isotropy of these systems ([Moore et al.l 
I1999T ). Our model suggests that the baryonic tracers 
are deeply embedded in their dark matter potentials to- 
day, resulting in higher values of Vn iax , proportional to 
2.2<7io S The difference is even more pronounced when 
comparing to the highest V m!iX values these systems ever 
achieved. Shown in the Figure as downward pointing 
triangles, these values indicate the amount of tidal mass 
loss that the subhalos hosting dSphs have undergone be- 
tween infall and today. In Figure we re-iterate this 
point by plotting the time evolution of Vmax as a func- 
tion of satellites' present-day luminosities. 

Finally, we show how VLI I dSphs compare to t he Tully- 
Fishcr relation (data from iPizagno et al.l I2007f ) for disk 
galaxies (Fig. 21 bottom right panel). Just like the MW 
dSph population, our star formation prescription fails to 
produce dwarfs that lie on the extrapolated Tully-Fishcr 
relation. The flattening of the distribution of spheroidal 
galaxies to a constant velocity value in comparison to 
field spiral galaxies is reminiscent of the "common mass 
scale" seen in the M300 mass function (see Fig. |9|), in 
which all currently detectable dSphs cluster around the 
value M300 ~ 10 7 Mq. As we shall see below, in our 
model the "common mass scale" is the result of a lumi- 
nosity bias, and future detection of fainter dSph satellites 
would break it and reveal lighter dark matter subhalos 
hosting stellar systems. 

We now turn to the structural properties of the surviv- 
ing VLII dwarfs. For all satellites observable in an all-sky 
SDSS-like survey, we calculate the circularly-averaged 
surface brightness seen by a galactocentric observer and 



4 This ~30% difference in Vm ax corresponds to a factor of 2 
increase in total mass, M <x V^* . 
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Fig. 3. — Left panel: Luminosity function of surviving satellites out to 280 kpc from the host halo center (shaded in green). Only galaxies 
that satisfy the magnitude-distance visibility SDSS criterion (see main text) are included. The shaded region for satellites fainter than 
My = — 11 spans the difference in number of dSphs falling within the SDSS DR5 footp rint from 2,000 rando mly oriented sky projections 
of the VLII halo Black dots: differential luminosity function of known MW satellites (sec Koposov ct al. 2009). Error bars arc indicative of 
Poisson number statistics in each magnitude bin. The lack of Magellanic Cloud equivalents in the VLII halo is apparent in the bright-end 
discrepancy between model and data. Note the different scales for the number of satellites fainter (left y-axis) and brighter (right y-axis) 
than My = —11, introduced to produce a continuous luminosity function despite the SDSS sky coverage of /drs = 0.194. Right panel: 
Evolution of the adopted mctallicity-luminosity relation after tidal stripping and dimming of the stellar populations with time. The original 
prescription (eq. [2{ assuming a mass- to-light ratio of unity is shown by a simple power-law (solid blue line). The stripping of stars from 
infalling satellites and the aging of stellar populations shifts luminosity to the left at a fixed value of [Fe/H]. Red diam onds: VLII satellites. 
The dashed and dotted lines show the obs erved metallicity-luminosity relation and la spread from lKirbv et ahl 1)20081 ). Black circles: Milky 
Way dSphs; data from Kirby et al. (2008) and references therein. 



fit it with a 2-D projected King profile (|Kindll962l ). 



h(R) = k 



^2 



-1/2 



- 1 



-1/2- 



(4) 



where is the surface brightness at radius R from 

the satellite center, k is the profile normalization, is 
the core radius of the satellite, and r t is its tidal (cutoff) 
radius. This formula provides a good fit to the surface 
brightness profiles of MW dSphs (|Mateol[l998l ). and has 
the advantage of having an analytical de-projected 3-D 
form in spherical geometry. Our use of the King pro- 
file as a fit to the projected surface brightness of VLII 
satellites instead of, say an exponential profile, is moti- 
vated by its standard use in published studies of dSph 
mass estimation. We find that the luminous component 
of the observable simulated dSph is well fit by a King 
profile across the full range of observable luminosities. 
Note that, while the tidal radii in the luminous compo- 
nent do not correspond to the physical truncation scale 
in the satellites' dark matter extent, the two scales are 
indeed comparable. On average, the ratio between the 
tidal radius in the luminous component and the tidal ra- 
dius in the dark matter extent varies between 0.5 and 
1 for satellites more luminous than 10 6 L Q , while it is 
between 0.3 and 0.8 for fainter systems. 

A comparison between the surface brightness and line- 
of-sight velocity dispersion profiles of two dSph satel- 
lites (Sculptor and Carina) and individual VLII coun- 
terparts is shown in Figure |6] Both surface brightness 
and kinematics should be compared only at radii larger 
than twice the VLII softening length of 40 pc, as the 
artificial "heating" of particles inside this region from 
numerical resolution effects would give a poor match to 



the data. In the surface brightness plot, we show the 
diverse light profiles of all VLII dwarfs that have a lumi- 
nosity within ±80% of the corresponding MW satellite. 
Equally diverse is the span of velocity dispersion profiles 
these dwarfs have, with values ranging between 4kms _1 
and 15kms _1 , pointing to the differences between the 
potential wells in which these systems reside. The in- 
sets in the figure show the evolution of the luminous sys- 
tems' parameters (King profile core and tidal radii and 
half-light radii) from 'tagging' (infall) until today. The 
behavior seen here is representative of the other SDSS 
detectable dwarfs in VLII - the core radii increase with 
time (though they remain close to the simulation soften- 
ing length), while the half-light and tidal radii decrease 
as the satellites lose stars from their outskirts to the main 
halo. In the velocity dispersion plots, we try to repro- 
duce the 'noisy' observed profile by sparsely sampling 
the simulation data. Using the same number of stellar 
particles in each corresponding bin as the one in the ob- 
served profile (M. Walker, private communication), we 
note that the scatter in the simulated profile increases 
significantly. 

In a systematic analysis, iStrigari et al.l (|2010| ) have 
shown that at least one simultaneously structural and 
kinematic counterpart to each of the five classical dSph 
with good kinematic and photometric data can be found 
in the Aquarius simulations. Our ability to obtain si- 
multaneous brightness and kinematic fits for dSphs like 
Sculptor and Carina can be seen as a confirmation of 
their results and a success of our model, and is a pre- 
requisite for testing the "common mass scale" of MW 
satellites. Structural parameters of the VLII satellite 
population (half-light radii and velocity dispersions mea- 
sured directly from the particle data, and central surface 
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Fig. 4. — Le/t <op panel: Cumulative distribution of galactocentric distance of VLII satellites fainter than My > —12. corresponding to 
the dSph satellites discovered by the SDSS. Black circles: Milky Way dSphs. Red diamonds: VLII satellites. Right top panel: Evolution of 
the stellar mass-halo mass relation from the infall epoch to the present in VLII surviving subhalos. Red diamonds: satellites detectable in an 
all-sky SDSS-like survey. Grey dots: all remaining (currently undetectable) satellites. The apparent increase in the star formation efficiency 
is caused by the preferential tidal stripping of the more diffuse dark matter component relative to the inner, more tightly bound stellar 
particles. Left bottom panel: Correspondence between 1-D line-of-sight velocity dispersion <ti os and 3-D maximum circular velocity V max 
for observable VLII dSphs. Solid blue line: Vmax tx V3a\ os . Red diamonds: VLII detectable satellites. Error bars indicate la variation 
over the 11 measurements of each profile used in our kinematic study. Dashed line: Vmax °c 2.2<ti os . Black triangles: Highest value of V ma x 
ever achieved by each satellite. Right bottom panel: The V ma x— I band magnitude relation for s urviving satellites (r ed diamonds) compared 
to to the Tully-Fisher relation of a broadly selected galaxy sample (black dots) from the SDSS (Pizagno ct al. 2007). The dashed line gives 
the extrapolation of the best fit to the galaxy data. VLII dSphs clearly do not follow the relation. 



brightnesses from the King model fits as a function of ab- 
solute magnitude) are compared to the observed proper- 
ties of MW dSphs in Figure [7] The agreement is reason- 
ably good. Reproducing the properties of the luminous 
component of dSph galaxies is a necessary condition for 
generating realistic "mock" data to test the mass mea- 
surement technique. 



3. KINEMATICS OF DWARF SPHEROIDAL SATELLITES 

3.1. The spherical Jeans equation 

The Jeans equation is the first moment over velocity of 
the collisionless Boltzmann equation. It provides a link 
between the underlying potential of a pressure-supported 
system and the distribution function of a tracer popula- 
tion in equilibrium with it. In spherical geometry its 



radial component is: 



dr 



(5) 



where v*(r) is the luminosity density profile of the 
tracer stars, ay is their radial velocity dispersion, and 
(3 = 1 — of /of i s tnc velocity anisotropy, i.e. the 
difference between the radial, ay, and the tangential, 
of = (<jg +ct|)/2, components of the velocity dispersion. 
Here, V c (r) is the circular velocity profile of the system, 
a measure of the underlying gravitational potential $, 



V c (r) 2 = G 



M(< r) 



dr 



(6) 



Since 3-D kinematic data for MW dSph stars are not 
available, one has to use a 1-D projection of the spher- 
ical Jeans equation along the line-of-sight between the 
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Fig. 5. — Comparison between the values of the 3-D maximum 
circular velocity V ma x at infall and today as a function of the satel- 
lite's present-day luminosity, for observable VLII dSphs. Red dia- 
monds: values of V m ax today. Black triangles: the highest values 
of Vmax ever achieved by each satellite. For the 10 most luminous 
satellites today, the V ma x values have decreased by a factor of 2.5 
on average (and up to 6 in individual cases) since accretion. 

observer and the studied system. This introduces an ad- 
ditional geometric effect inv olving the velocity anisotropy 
(jBinnev fc Tremaind l200l : 



' las 



S (R) = 



v^a^r 



s/r 



R 2 



dr, (7) 



where is the projected surface brightness profile 

of the tracer population. Of all quantities involved in 
the projected Jeans equation, a\ os and I*(R) are readily 
determined observationally, and v*(r) can be obtained 
through a dc-projection of the surface brightness profile, 
ass uming i t follo ws so me analyti c al form - either a spher- 
ical Kind (|1962f > or a iPlummerl (|1911l ) profile. Both of 
these profiles adequately fit the light distribution of MW 
satellites and the choice of one or th e other does not affec t 
the results of our kinematic study (jStrigari et al J 120081 ) . 
The remaining parameters, namely /3(r) and V c {r), are 
not constrained observationally and have to be modeled 
appropriately. 

We adopt a generalized NFW profile (INavarro et al.l 
fl996Tl for the mass distribution of subhalos: 



p(r) 



(r/r,)T(l + r/r s ) 



s ' 



(8) 



where p s gives the density normalization, r s is the scale 
radius, and 7 and S define the inner and outer slopes of 
the density fall-off. The fiducial NFW form has 7=1 
and 5 = 2. To allow for a range of inner density cusps 
as well as outer tidally-modified profiles, we allow these 
two parameters to vary in the inter vals 0.7 < 7 < 1.2 
and 2 < S < 3 (jStrigari et all 120071) . We have verified 
that the upper limit on the outer slope S is high enough 
to accommodate the profiles of even the most massive 
VLII subhalos, whose density profiles fall off faster than 
those of field halos in the tidally-stripped outer regions 
(|Diemand et al.ll2008l) . 

The velocity anisotropy parameter (3 is largely uncon- 
strained by current observations. According to its defini- 
tion, /3 can vary between —00 and 1, with /3 = —00 when 



Ctot = ot an< i P = 1 when er t ot = <?>■ Available data are 
often in the form of 1-D line-of-sight velocity dispersion, 
insufficient to constrain ay and at individually. These 
two quantities are well determined only at the very cen- 
ter and very edge of a system, since a\ os (0) = a r and 
cios(^t) — <?t- In between these two extremes, the mass 
of a dSph enclosed within a radius close to the system's 
half-light radius can be calculated such that the uncer- 
tainty from the velocity anisotropy term is minimized 
(jWolf et al.H20Tl . At any other intermediate radius, the 
(3 parameter will more strongly affect the mass measure- 
ment and must be estimated using a theoretical model. 
In principle, /3 can be a function o f radius within a partic- 
ular system (jStrigari et al.l l2007i). Previous studies have 
shown that marginalizing over a set of parameters that 
define a general anisotropy profile does not constrain the 
mass distrib ution better than as suming (3 to be constant 
with radius (jWalker et al.l [20 09). The latter case, how- 
ever, grea tly simplifies the solutio n to the spherical Jeans 
equation (jMamon & Lokas 2005): 

p OO 

v*a 2 r (R) = GR~ 2p r 2 ^ 2 ^(r)M(r)dr. (9) 



The standard procedure in mass modeling of pressure- 
supported systems is to generate a theoretical veloc- 
ity dispersion profile that fits the observed one, and to 
marginalize the solution over all free parameters but the 
density profile normalization. While most model param- 
eters are largely unconstrained by this fitting procedure, 
the density profile normalization (defined as the total 
amount of mass enclosed in a fixed aperture radius, typ- 
ically 300 pc or 600 pc) i s very well determi ned by the 
data (jStrigari et al.l 120081: IWalker et all 12009( 1 . It is this 
analysis that yields the so-called "common mass scale" 
of MW dSph satellites. 

3.2. Applications to VLII dwarf satellites 

We apply the mass measurement technique above to 
the simulated data using the following procedure: 

1. Obtain a\ os (R) and I*(R) for a satellite of interest; 

2. Choose the density profile and anisotropy parame- 
ters p s , r s , 7, 8, and j3; 

3. Integrate the Jeans equation (JSJ) to solve for ay(r); 

4. Project a r (r) along the line of sight (eq. [7]); 

5. Compare the solution to the measured a\ os {R)] 

6. Iterate (2 - 5) to find the model that best matches 
the kinematic data. 

Once the best fit mass density profile is determined, 
we integrate it to a given aperture radius, either 300 pc 
or 600 pc from the satellite center to obtain the total 
enclosed mass, M300 or Mgoo, and compare it directly to 
the "true" value in the simulation. To obtain the "mock" 
observations of VLII dSphs in step 1, we place "hypo- 
thetical observers" in the simulation volume, and reduce 
computational time by "observing" each of the 65 de- 
tectable satellites only 11 times. One of the "observers" 
is always in the center of the VLII main host, while the 
other 10 are randomly distributed at the same distance 
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Fig. 6. — Surface brightness and velocity dispersion profiles of two VLII satellites (red dots) closely resembling the Milky Way dSph 
Sculptor and Carina. Error bars are indicative of la Bootstrap re-sampling errors in each bin. Simulated surface brightness and kinematics 
should be compared only for R > 80 pc, corresponding to twice the VLII softening length. The observed surface brightness prof iles are 
plotte d with dashed lines and the uncertainty is shown by the shaded region. Brightness fits for Sculptor and Carina are taken from lMateol 
1119980 . while the kinematic data (black dots with error bars) arc from Walker ct al. (2009). Also shown is the diversity in brightness profiles 
among VLII satellites with luminosities within ±80% of Sculptor's and Carina's (grey dot-dashed lines). The insets follow the evolution of 
the King profile structural parameters for the fiducial satellites from the time of 'tagging' (infall) until today. Solid blue lines: the 95% 
confidence region in the kinematic profiles where the data are sparsely sampled as the observations. 



from the satellite as the galactocentric "observer". In 
this way we attempt to capture the anisotropy of subha- 
los' physical shapes and velocity ellipsoids. 

Since integrating the Jeans equation requires a bound- 
ary condition on a r (r t ), which is not known a priori from 
the data, the line-of-sight projected solution cannot be 
trusted across the full extent of the surface brightness 
fit. The observed line-of-sight velocity dispersion profile 
is therefore calculated in 10 linear distance bins out to 
60% of each satellite's tidal radius r t . The simulated data 
used to determine each subhalo's kinematic profile do not 
have an associated observational uncertainty. To make 
a fair comparison with the derived mass estimates, "ob- 
servational" errors sho uld be consistent w ith the ones in 
the current literature (|Walker et al.ll2009t) . This is done 
by randomly assigning error values <ii ranging between 1 
and 2 kms -1 to the velocity dispersion data when eval- 
uating the goodness of fit in equation (jTTJ|) below. 

Sample line-of-sight velocity dispersion profiles for two 
VLII simulated satellites and their assigned "observa- 
tional" errors are shown in the top left panel of Figure [5J 
The two satellites have masses of M 300 = 7.22 x 10 6 M Q 
and M 30 o = 2.65 x 10 7 M Q , respectively. For each satel- 



lite, we plot all 11 best fit cri os (i?) profiles. The two sets 
of data points (in black and grey) correspond to the kine- 
matic profile fits highlighted in darker colors, which were 
chosen to the span the range of M300 for which the fit's 
likelihood falls to 10% of its peak value (see the right 
top panel of Fig. [5]). The average M 300 values inferred 
from these 11 fits are 6.65 x 10 6 M and 2.44 x 10 7 M Q 
respectively, within 8% of the true values. 

To get the above estimates, we fit a model to each "ob- 
served" kinematic profile using a Monte Carlo Markov 
Chain (MCMC). A theoretical velocity dispersion profile 
cios,th(-R) is computed for each parameter set in step 2. 
and the likelihood of the fit is evaluated as: 



1 



v/2^ 



:CXp 



1 (a\ os (Ri) - erios,th(-Ri)) 5 



(10) 



where <ri os (Ri) is the line-of-sight velocity dispersion 
measured at radius Ri and <Ji is its associated error. Each 
MCMC chain is initialized to maximize the likelihood C 
on a coarse model grid that spans the full parameter 
space. It is then left to explore the parameter space fol- 
lowing the Metropolis-Hastings algorithm. Two concur- 
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Fig. 7. — Structural properties of surviving VLII satellites. Plotted against V-band absolute magnitude are their half-light radii (left 
panel), line-of-sight velocity dispersions (middle panel), and central surface brightnesses (right panel). Red diamonds : observable VLII 
dwarfs. Black circles with error bars: Mil ky Way dSph. Data compiled as follows: My. classical dw arfs: (Mateo 1998), except for L eo II 
llColeman et al.H20q7l ) ultra -faint dwarfs: UMartin et'a l. 2008), except for Hercules ISand et alJI2009Q an d Leo T flde Jong ct al. 200 8)); r e g 
and p-i^T T^oTfet*^TTll2010l) : /ioy: classical dwarfs; {Mateo 199?|), except for Draco lMartTri^tani2*008l ) , ultra- faint dwarfs: IIMartin et aT] 
I200j ), except for Leo T {Ir win et al. 2007). The dash-dot red vertical line in the left panel marks 2.8x the VLII softening length of 40 pc 
fsee lCooper et al.ll2010l) . 



rent chains are run on each data set and the MCMC is 
stopped after 10,000 attempts of each chain. This num- 
ber is large eno ugh to provide sat isfactory fits to the 
kinematic data (| Walker et al.1 [20091 ) . The 11 resulting 
likelihood functions for M 30 o of the two satellites above 
are shown in the right top panel in Figure [8] 

Once the best-fit model for each data set has been 
found, it provides the full mass density profile that can 
be integrated to any distance from the dSph center, with 
300 pc or 600 pc being of particular interest. There are 
11 likelihood distributions for the estimates of M300 and 
Mgoo of each VLII satellite, and these can be compared 
to the true values. The results are shown in the bottom 
panels of Figure [8j in which wc plot the average of the 
11 best mass estimates as a function of the true dark 
matter mass. The nature of the MCMC data allows us 
to make two estimates for the error in each case. The 
smaller error bars (in black) in Figure [8] reflect the span 
of the 11 best fit values for each satellite. The larger 
errors (in grey) span the full range in values where the 
likelihood of the fit has fallen to 10% of its maximum 
value. For the particular two satellites in the top panels 
of Figure [51 wc thus have the following mass estimates: 
M 300 = 6.65±g;|4±|;|g x 10 6 M Q for the lower-mass satel- 
lite and M300 = 2-44±o:il±\ir x 10 7 Af Q for the higher- 
mass one, where the first set of errors corresponds to the 
best fit value range and the second set of errors corre- 
sponds to the 10% likelihood value range. 

We infer that the spherical Jeans analysis is a robust 
method to estimate the central masses of dSph satel- 
lite galaxies. The mean best fit A/300 values determined 
through the Jeans analysis are on average within 12% 
of the true VLII values. The normalized scatter among 
the 11 best fits about the mean value is 15 on aver- 
age and is shown as solid lines in the bottom panels of 
Figure [5] The normalized full range scatter about the 
mean where the likelihood of the fit is greater than 10% 
of its peak value is io'53 on average (dashed lines). No 



significant difference to these results arises when just the 
galactocentric observer's estimates are considered: the 
best fit galactocentric estimate of M300 across the full 
mass range is on average within 16% of the true VLII 
value. 

When evaluating M 6 oo, we exclude systems that have 
Mqqo < 1.5 x 10 7 M . It is apparent in the bottom right 
panel of Figure that these exhibit a larger scatter than 
the systems with Afgoo > 1-5 x 10 7 M Q . The smallest 
dSphs host stellar systems smaller than 0.6 kpc in radial 
extent. The estimate of Mqqq for them is therefore an ex- 
trapolation to a region in parameter space where tracers 
are simply not present. Among the dSphs with kinematic 
data available all the way to r = 0.6 kpc, the mean best 
fit Afgoo value determined with the Jeans analysis is on 
average within 17% of the true VLII value. The normal- 
ized scatter among the 11 best fits about the mean value 
is io'26 on average. The normalized full range scatter 
about the mean where the likelihood of the fit is greater 
than 10% of its peak value is 55 on average. The VLII 
galactocentric observer would report a A/goo value within 
23% of the true VLII value. 

The range in ai os (R) and 7*(i?) when measured from 
different viewpoints is due to the anisotropic velocity ten- 
sor of the dSph and any deviation from perfect spheric- 
ity of its shape. The mass estimate uncertainties among 
the best-fit values presented here should therefore be re- 
garded as intrinsic, cosmologically-induced ones that will 
eventually be minimized when 3-D kinematic data is col- 
lected and better models of the nature of dSph satel- 
lites in high-resolution, fully cosmological hydrodynami- 
cal simulations are carried out. 

In the left panel of Figure we plot the mass en- 
closed in the inner 300 pc versus total luminosity for all 
satellites that would be visible according to the SPSS 
detect ion limit. Ovcrplottcd arc data from Strigari et al. 
([2008ll . suggesting a "common mass scale" . Observable 
VLII dSphs exhibit a behavior similar to the MW dSph 
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Fig. 8. — Left top panel: Sample line-of-sight velocity dispersion profiles for two VLII dSph, of masses M300 = 7.22 X 10 6 M Q and 
M300 = 2.65 x 10 7 M . Best fit Jeans solutions for cri os (R) plotted in red and blue. Fits highlighted in darker colors correspond to the 
points with error bars and the envelope 10% likelihood functions shown in the right top panel. The average best fit estimates of M300 in 
this case are 6.65 X 10 6 Mq and 2.44 X 10 7 Mq respectively. Right top panel: Likelihood functions of M300 for the two satellites in the 
left panel. The 11 realizations span a range in the best fit values around their mean (black dots with error bars) as well as a range in 
uncertainty where the likelihood has fallen to 10% of the best fit value (grey error bars). The vertical dot-dash lines indicate the true VLII 
values of JV/3oo- Left bottom panel: Total mass within the inner 300 pc, M300 as determined through the Jeans equation in comparison 
with its true value from VLII. Mean values of 11 individual measurements arc plotted with black points. Black error bars span the range 
between best fit values, while the grey error bars span the full range of values where the likelihood of the fit has fallen to 10% of its peak 
value. Solid red line: Unity correspondence. Solid blue lines: Average over the full M300 mass range of the spread among the best-fit 
measurements. Dotted blue lines: Average over the full M300 mass range of the spread among measurements where the likelihood is 10% 
of its peak. Right bottom panel: Total mass within the inner 600 pc, M^oo as determined through the Jeans equation in comparison with 
its true value from VLII. Labels as in the M300 panel on the left. Systems with A/600 ^ 1.5 x 10 7 Mq (left of the red dash-dot line) do not 
have kinematic data out to r = 600 pc, which results in a larger scatter in their estimated Me,oo values. 



population, though a linear fit suggests a shallow but sig- 
nificant correlation between M300 and the total luminos- 
ity: M300 oc ^o ossiO- 024 ^ where the error on the slope re- 
flects the la deviation from bootstrap re-sampling of the 
sample. This is steeper than the £°- 03 ±°- 03 relat ion de- 
rived from observational data (jStrigari et a l. 2008). Low- 
ering the derived M^nn values of some of the ultra-faint 
dSphs, as suggested bv lAden et al.l (|2009f ) for the case of 
the Hercules dwarf, might be sufficient to reconcile the 
observations with the slope predicted from our model. 

In our model, the nearly constant M300 is explained 
by the combination of the steep adopted stellar mass- 
halo mass relation (eq. [T]) and the SDSS detectability 
magnitude limit. The remaining ~ 2, 000 VLII surviving 



subhalos to which we assigned stars (about 200 of which 
are more luminous than My = 0) are just too faint to be 
detected by a SDSS-likc survey. Some of them have M 300 
significantly lower than 10 7 Mq (see the right panel of 

F igure H). 

iBovlan-Kolchin et al.l (|201 If ) has recently pointed out 
that state-of-the-art collisionless simulations of Galactic 
structure, both VLII and Aquarius, predict a population 
of massive subhalos that appear to be too dense (i.e. 
have too high a Vmax for a given r vmax ) to host any of 
the known MW dSph. While a failure of the CDM model 
is one interpretation of this result, another possibility is 
that star formation becomes stochastic below some mass, 
in which case the MW halo would harbor several massive 
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L [L ] M 300 [M ] 

Fig. 9. — Left panel: Total mass within the central 300 pc, M300 as a function of total luminosity for dwarf satellites in VLII visible by 
an all-sky SDSS-like survey (black dots). Currently undetectable VLII d warfs are shown in gr ey. A linear best fit to the data with slope 
0.088 is shown with a red solid line. Data points with error bars are from:Strigari ct al. (2008), indicating a "common mass scale". These 
error bars correspond to the points where the likelihood is 60.6% of its peak value on either side of the best estimate. For comparison, we 
plot our average error estimates in M300 from Figure \8\ about the best fit (solid red line): solid blue lines: full spread in best-fit estimates; 
dotted blue lines: full spread in values where the likelihood is 10% of its peak. The "classical" dSphs (L > 2 X 10 5 Lp^) are shown as blue 
squares and the "ultra-faint" dSphs as red circles. The improv ed mass for the Hercu l es dw arf from [Aden et a l. (2009) is shown as a red 
upward pointing triangle. The six VLII subhalos identified by Boylan-Kolchin et al. (2011) as too dense to host any of the known MW 
dSphs are plotted in orange as five circles (if detectable) and a single diamond (currently undetect able). Right panel: M300 mass function 
for VLII dSph satellites. Red diamonds: MW dSphs with masses taken from jStrigari ct al. 20Q3). Black dots: VLII satellites detectable 
in SDSS. Grey dots: all VLII surviving subhalos that were tagged with stars, including those currently undetectable. The "common mass 
scale" is apparent in the currently detectable population. 



but dark satellites. It is also quite plausible that baryonic 
feedback processes in the more massive subhalos redis- 
tribute some of their da rk matter, thereby lowe ring their 
central density (see e.g. iGovernato et al.1120101 ). For ref- 
erence we have marked th e six VLII subhalos that are 
too dense according to the iBovlan-Kolchin et al.l (|2011[ ) 
analysis in the left panel of Figure [9l Of these six, five 
have a luminosity greater than My = —11 and would 
have been detected as luminous dSph satellites. In our 
model, it is thus not possible for these overdense sub- 
halos to have evaded detection because of an effectively 
low star formation efficiency resulting, for example, from 
their accretion or orbital history. As expected, these six 
halos have the highest values of M 300 for their luminos- 
ity among all VLII satellites. A reduction in M300 by 
about a factor of two would bring them into agreement 
with the classical MW dSph sample. Our study of the 
standard Jeans analysis demonstrates that the 'observed' 
M300 values cannot be changed by that amount, i.e. that 
errors in the mass estimates are not responsible for the 
discrepancy between the models and the observations. 
Future, high-resolution hydrodynamical simulations of 
MW analogs that include an accurate treatment of bary- 
onic processes and their effects on the central densities of 
luminous satellites will be crucial to the settling of this 
outstanding issue. 

4. SUMMARY 

We have used a particle tagging technique to dynami- 
cally populate the iV-body Via Lactea II high- resolution 
simulation with stars. The method is calibrated us- 
ing the observed luminosity function of MW satellites, 
the concentration of their stellar populations, and the 
metallicity-luminosity correlation they exhibit, and sclf- 
consistently follows the accretion and disruption of pro- 



genitor dwarfs and the build-up of the stellar halo. The 
particle tagging technique makes full use of 27 simula- 
tion snapshots spanning the assembly history of the VLII 
main host between z = 27.54 and z = 0. A subset of the 
most bound dark matter particles in each of about 3,200 
progenitor subhalos is tagged as stars at infall, when star 
particles acquire a stellar mass, an age, and a metallic- 
ity. Their subsequent evolution is purely photometric 
and kincmatical in character, as the collisionless stellar 
populations age following population synthesis models, 
and are accreted and disrupted in a comological "live 
host" . Luminous particles linked to self-bound substruc- 
tures in the z = snapshot constitute a surviving pop- 
ulation of dSph satellites. Lik e the previous studie s of 
IBullock fc Johnston! (|2005D and lCooper et al.l (|2010t) . we 
have shown that simple prescriptions for assigning stel- 
lar populations to collisionless particles are able to re- 
produce many properties of the MW dwarf satellite pop- 
ulation, like velocity dispersions, sizes, brightness pro- 
files, metallicities, and sp atial distribution. In agreement 
with previous results (e.g. |Kravtsov 20101 : IKoposov et al.l 
120091 : [Madau et al. 20083), a steep stellar mass-subhalo 
mass relation is required to fit the observed luminosity 
function of MW dSphs and offer an empirical solution 
to the missing satellite problem. This model predicts 
the existence of approximately 1,850 subhalos harbor- 
ing "extremely faint" satellites (with mass-to-light ratios 
> 5 x 10 3 ) lying beyond the SDSS detection threshold. 
Of these, about 20 are "first galaxies" , i.e. they formed 
a stellar mass above 10 M before redshift 9. The ten 
most luminous satellites (L > 10 6 L©) in the simulation 
are hosted by subhalos with peak circular velocities to- 
day in the range V max = 10 — 40 kms -1 that have shed 
between 80% and 99% of their dark mass after being ac- 
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creted at rcdshifts 1.7 < z < 4.6. The satellite maximum 
circular velocity Vi nax and stellar line-of-sight velocity 
dispersion cri os today follow the relation Vi nax = 2.2cri os , 
i.e. stellar particles in satellites have smaller velocity 
dispersions than the dark matter. The ability to obtain 
simultaneous brightness and kinematic fits for dSphs like 
Sculptor and Carina is a confirmation of earlier searches 
for MW satel lite counterpa r ts in simulated halos such 
as the one by iStrigari et al.l (|2010f ) and a success of our 
model. It is a prerequisite for studying another small- 
scale puzzle of galaxy formation in ACDM, namely the 
observed "common mass scale" of MW satellites, i.e. the 
fact that, despite spanning about 5 decades in luminos- 
ity, MW dSphs appear to contain a similar amount of 
dark mass, M 30 o ~ 10 7 M Q , within a fixed aperture ra- 
dius of 300 pc fro m their centers (jStrigari et al.1 [20081 : 
IWalker et al.ll200"9| ) . We have applied a standard mass es- 
timation algorithm based on Jeans modelling of the line- 
of-sight velocity dispersion profiles to the simulated VLII 
dwarf spheroidals, and tested the accuracy of this tech- 
nique. We have shown that the mass enclosed within the 
inner 300 pc can be estimated accurately to within 20% 
of the true value. The predicted Af3oo-luminosity rela- 
tion for currently detectable simulated satellites is nearly 
flat. We do, however, predict a weak, but significant pos- 
itive correlation for these objects: M300 oc £ 088 ± - 024 . 
Our results are consistent with previous studies in which 
satellites observable by presently available wide-field sur- 
veys naturally tend to reside in subhalos that contain 
about 10 7 Mq within their inner 300 pc. Future, even 
deeper surveys with the Large Synoptic Survey Telescope 
may be able to determine if an extension of this satel- 
lite galaxy population, at an even fainter magnitudes, 
indeed exists in the halo of the Milky Way. These exper- 



iments will cast light on the processes that govern galaxy 
formation on even smaller scales predicted in ACDM 
cosmology. Ultimately, the prescriptions used in this 
work must be tested with hydrodynamical simulations 
that include baryonic physics. Simulations of the for- 
mation of late-type spiral galaxies in a cold dark matter 
(ACDM) universe have traditionally failed to yield re- 
alistic candidates. Our group has recently reported a 
new cosmological hydrodynamic simulation of extreme 
dynamic range, "Eris" , in which a close analog of a 
Milky Way disk galaxy a rises naturally at the present 
epoch (jGuedes et al.ll20ll . With a mass resolution that 
is ~ 20 times lower than VLII, Eris includes radiative 
cooling, heating from a cosmic UV field and supernova 
explosions (blastwave feedback), a star formation recipe 
based on a high gas density threshold, and neglects any 
feedback from an active galactic nucleus. The study of 
the high-mass end of Eris' satellite population and the 
comparison with the results presented here will be the 
subject of a subsequent paper. 
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